/****************************************************************************************************
* Title: Workhorses of Opportunity, Howard and Weinstein
*
* Inputs:
*   - sf12000countydistance500miles.dta
*   - justnormasylum.dta
*
* Outputs:
*   - asylumdistance.dta
*   - closeasylums.dta
*   - close_far_asylums.pdf
*
* Description:
*   This script computes the distance from each asylum-only county to the nearest normal school county,
*   classifies "close" asylums using the median distance, generates a state-level summary, and maps the results.
****************************************************************************************************/

* Load county-to-county distance data and clean
clear
use "sf12000countydistance500miles"
destring county1, replace
destring county2, replace
tempfile distance
save `distance'

* Load normal school counties
clear
use "justnormasylum"
rename cty_fips county2
keep if hasnormalschool
tempfile normal
save `normal'

* Load asylum-only counties
clear
use "justnormasylum"
keep if !hasnormalschool
rename cty_fips county1
tempfile asylum
save `asylum'

* Merge distances from asylum counties to all normal school counties
merge 1:m county1 using `distance', keep(1 3) nogen
merge m:1 county2 using `normal', keep(2 3) nogen

* Find minimum distance to a normal school for each asylum county
collapse (min) mindistance = mi_to_county, by(county1)
rename county1 cty_fips
save "asylumdistance.dta", replace

* Use the median minimum distance as the cutoff for "close asylum"
summ mindistance, detail
gen closeasylum = mindistance < r(p50)

* Merge back with full county-level indicators
merge 1:1 cty_fips using "justnormasylum"
replace mindistance = 0 if hasnormalschool == 1

* Plot map of minimum distance to normal school
gen county = cty_fips
maptile mindistance, geo(county2000) mapif(mindistance != .) ///
    stateoutline(.1) cutv(0 41.05) fc(Spectral)
graph export "close_far_asylums.pdf", as(pdf) replace

* Save classification for later use
keep cty_fips closeasylum hasnormalorasylum hasnormalschool
save "closeasylums.dta", replace

* Tabulate states with variation across county types
replace closeasylum = . if hasnormalschool
gen state = floor(cty_fips / 1000)
gen farasylum = 1 - closeasylum
collapse (max) hasnormalschool farasylum closeasylum, by(state)
list

* Summary of states by presence of each type
summ state if hasnormalschool & farasylum == 1 & closeasylum == 1
summ state if hasnormalschool & farasylum == 1 & closeasylum == 0
summ state if hasnormalschool & closeasylum == 1 & farasylum == 0
